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ABSTRACT 

Cross-power spectrum is a quadratic estimator between two maps that can 
provide unbiased estimate of the underlying power spectrum of the correlated 
signals, which is therefore used for extracting the power spectrum in the WMAP 
data. In this paper we discuss the limit of cross-power spectrum and derive the 
residual from uncorrelated signal, which is the source of error in power spectrum 
extraction. We employ the estimator to extract window functions by crossing 
pairs of extragalactic point sources. We desmonstrate its usefulness in WMAP 
Difference Assembly maps where the window functions are measured via Jupiter 
and then extract the window functions of the 5 WMAP frequency band maps. 

Subject headings: cosmology: cosmic microwave background — cosmology: ob- 
servations — methods: data analysis 

1. Introduction 

Variables and parameters in theories of cosmology are often derived based on a statis- 
tical ensemble. In practice, particularly in data analysis on the CMB, however, one often 
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encounters the well known and rather unique issue, namely we can do statistical analysis 
on only one version of the celestial sphere using an ergodic hypothesis. For example, the 
angular power spectrum of the temperature anisotropies Ci = (a£ m a} m ), but with one version 
of the sphere it can be estimated with only 2£ + 1 modes: C/ St = (2£ + Ylm=-e a im a \m- 
Lack of a statistical ensemble thus induces in power spectrum analysis the uncertainty called 
"cosmic variance" (White, Krauss and Silk 1993). Another issue arising from lack of an en- 
semble is coined "cosmic covariance" , which refers to a non-zero chance correlation between 
the CMB and the foregrounds (Chiang, Naselsky and Coles 2009). The cosmic covariance 
is the source of error in employing quadratic minimization for separating the CMB and the 
foregrounds (Bennett et al. 2003; Hinshaw et al. 2007). 

In this paper we address another issue in cross-power spectrum arising from lack of 
an ensemble. Cross-power spectrum is useful in eliminating uncorrelated signal and is thus 
commonly used in power spectrum estimation from multifrequency measurement of CMB. It 
is required, however, to have an ensemble of modes in the summation to completely eliminate 
the uncorrelated signal. Therefore lack of an ensemble of modes results in residual (hence 
error) in the cross-power spectrum. 

Another important issues related to extraction of the power spectrum of CMB tempera- 
ture anisotrpies is the calibration of window function convolving the signal, which is usually 
by measuring Jupiter (Page et al. 2003). It can also be via measuring bright extragalactic 
point sources as they are manifestation of the window function as well, albeit with back- 
ground CMB and instrument noise. We can then employ cross-power spectrum on pairs of 
extragalactic point sources from different patches of the sky, where the background CMB 
and other signals are uncorrelated and can be largely eliminated. 

This paper is arranged as follows. In Section 2 we examine the limit of cross-power 
spectrum as a quadratic estimator. We in Section 3 employ the estimator on a pair of 
simulated point sources to obtain the window function. In Section 4 we put this method 
to the test on WMAP Difference Assembly (DA) maps whose window functions are already 
carefully calibrated via measuring Jupiter and we extract the window functions of WMAP 
frequency band maps. The conclusion and discussion is in Section 5. 
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2. Cross-power spectrum and its limit 

Cross-power spectrum (XPS) is a quadratic estimator between two maps j and f, and 
on a sphere it is written as 

1 l 

m=— £ 

where je m and j' £m are the spherical harmonic coefficients of the two maps, and * denotes 
complex conjugate (Hinshaw et al. 2003). In order for the XPS to be real, the usual practice 
is by replacing the j\ m j\ m in Eq.(l) with 

titrJtm + j'Ljtm)/ 2 = \jirn | \j' tm | COS A0 £m , (2) 

where A0^ m is the phase difference between j £rn and j' £m . The advantage of XPS as an 
unbiased quadratic estimator for power spectrum estimation lies in the fact that XPS returns 
with its usual power spectrum \je m \ 2 if j and j' are of the same signal, and if, on the other 
hand, j and j' are uncorrected then A(fie m is distributed uniformly random in [0, 27r] and 
thus (%2 m= _ e \3im\\j'trn\ cos A^fot) = 0, where the angle brackets denote ensemble average. 
Hence XPS is useful in eliminating uncorrelated signals while preserving the correlated one 
and has been employed by WMAP to extract CMB spectrum by crossing the Differencing 
Assemblies (DA) (Hinshaw et al. 2003, 2007; Nolta et al. 2009; Larson et al. 2010) . 
In practice, however, one does not have a statistical ensemble and the capability of XPS 
to extract the underlying power spectrum from correlated signal depends on how much the 
uncorrelated signal is decreased with finite realizations and finite number of multipole modes 
involved in the summation at £. 

Take as an example extracting the CMB signal from WMAP data. The Wl and W2 
DA have the same CMB but uncorrelated noise, so writing a™ = a\ m + n^ m and a^ 2 = 
a C £ m +n' em , in XPS the correlated signal \a c lm \ 2 is the desired one whereas those uncorrelated 
terms between CMB and noises AV n , X e CQ and between noises X/ n shall be decreased but 
nevertheless not to zero, which is then the source of error in the XPS. 

In order to examine the level of residual in XPS, we consider two Gaussian and uncor- 
related maps j and j', whose spherical harmonic coefficient and j' em have zero mean and 
variance and .][ respectively and the corresponding phase difference A(f>e m is uniformly 
random in [0,27r]. One should note that the distribution of Eq.(2) under this condition is a 
normal product distribution, which is itself not Gaussian but can be expressed as 
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where K is the zeroth-ordered modified Bessel function of the 2nd kind (Springer 1979). 
Although the distribution of the variable is not Gaussian, what we need to examine is its 
summation. Since the XPS defined in Eq.(2) can be either positive or negative, we look 



instead at its root mean square y ((X/ J ') 2 ). Due to equal partiution on both perpendicular 
directions of complex plane, it can be written: 



(4) 



the right hand side of which is closely related to a two-dimensional random walk. Random 
walk statistics is employed as a Gaussianity test on CMB (Naselsky et al. 2005; Stannard 
& Coles 2005). Recall that the celebrated Pearson's random walk is a 2D isotropic random 
walk with equal step-length, expressed in the complex plane r = ^ k \dk\ exp(i<pk) where 
\dk\ — d. The pdf of the resultant displacement r = \r\ after £ steps has the form 



p(r) = ^ ex H~M^' (5) 

from which (r 2 ) = £d 2 . Pearson's walk can be generalized to unequal step-length and it can 
be shown (Pearson 1905; Rayleigh 1905) that when £ is large and the mean-square step 
length is known : (f) -1 \dk\ 2 = d 2 , the pdf for the displacement has the same form as 
Eq.(5) (Hughes 1995). For unequal step length \d k \ = \jem\\je m \, the ensembled-average 
square step length is then (|j^m| 2 |j^ m | 2 ) = JtJ[i so the mean square displacement of the 
random walk (| Y m =-e \jem\\je m \ exp(iA0^ m )| 2 ) = (21 + l)J^J' t Note that in driving the 
result we only assume j and f are Gaussian, so they can be either white noise (with flat 
power spectrum) or others whose power spectrum has scale dependence. In any case, we can 
write 

<(X^') 2 > ! 

(6) 



yfjj* 2y/e+l/2 

The decreasing of the "junk" is inversely proportional to the square root of the number of 
random walk, and it can be further decreased by l/y/LN with binning L = A£ multipole 
numbers and averaging from iV sets of XPS. We demonstrate Eq.(6) in Fig.l by plotting 
Xi from two white-noise maps. One can see the the decreasing of the power spectrum 
proportional to (2^Jt + 1/2) -1 . Note that Xe is jumping either postively or negatively and 
only positive points are plotted. 

Knowing the residual of the junk from XPS, we can estimate its capability in extracting 
CMB. In crossing iV pairs of maps containing CMB with the same noise level, we assume 
the CMB power spectrum is C e ~ A£~ 2 and the instrument noise level is N e = RA where 



1 10 100 1000 

I 



Fig. 1. — Demonstration of Eq.(6). The black dots denote XPS from two white noise maps 
with unity power spectrum, whereas the larger dots are the average from 200 pairs of white 
noise maps, binned with L = Al = 50. The upper line is (2^/£ + 1/2) -1 and the lower 
{2y/LN(£ + 1/2)) -1 , where L = 50 and iV = 200. Note that the negative values are not 
plotted. 
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R indicates the noise level relative to the CMB and R <C I 1 . The residual level from 
between CMB and noise, and between noises are 2X/ n ~ Ay/R(£ 3 LN)~ 1 ^ 2 and X/ m ' ~ 
RA{ULN)~ 1 / 2 , respectively. They intersect at £ ~ 2/^/R. For £ < 2/VR, X e cn is the main 
residual with contribution: 2X/ n /C/ ~ a/ Ri/LN, and LA > £ 2 i? in order for the residual 
contribution to be within cosmic variance, which shall be the case as long as LN > 4. For 
£ > 2/\/R, on the other hand, the residual in the CMB extraction is R^J £ 3 / ALN , and 
beyond £ = (ALN/R 2 ) 1 ^ the residual level is higher than the CMB, thereby XPS fails. 

For a more realistic situation, where the XPS is employed on signals that is convolved 
with window function: A£~ 2 exp(— B 2 £ 2 ), we plot in Fig.2 the upper limit of R (the noise 
level) if the XPS is to retrieve the CMB at the corresponding multipole number. The different 
FWHM are chosen according to those from channels of ESA Planck Surveyor. As expected, 
more binning and more sets of XPS shall retrieve higher modes for the same noise level. 

In Figure 3 we show the XPS (dots) from two simulated maps. Top panel shows Gaussian 
random signal with power spectrum Ce = 10 4 £~ 2 (/iK 2 ) convolved with beam FWHM 30 
arcmin in both maps, the theoretical power spectrum of which is shown in dash-dot curve. 
They are added with pixel (white) noise signals with differnt levels N e = 1CT 1 and N' e = 
10" 3 (/iK 2 ). The residual from XPS (i.e. Eq.(6)) is plotted in dash line. Note that the 
residual in employing XPS on maps with different noise levels is proportional to a/ N^N' e , 
which reflects the real situation of crossing between different frequency bands such as WMAP 
V and W band. In bottom panel we simulate 4 WMAP W channel DA maps, convolving 
with beam FWHM 12 arcmin before adding WMAP simulated noise. There are 6 XPS and 
the binned (A£ = 10) retrieved power spectrum is denoted with plus sign. One can see the 
residual from the noise is the main error from XPS, even if it's not white noise. 

3. Application of cross power spectrum 

One of the applications of XPS is to retrieve window function from bright extragalactic 
point sources. The standard way of measuring the window function in CMB experiments is 
via measuring planets such as Jupiter (Page et al. 2003; Hill et al. 2009). 

The importance of window function can be illustrated in the following simple example. 
The final step of CMB power spectrum estimate often involves deconvolution as the retrieved 
singal (either via ILC method or foreground template fitting) is also convolved with a com- 



1 The noise level Nt = RA is related to noise fluctuation cr 2 through a 2 — N p i x RA/4Tr, where iVpj x is the 
total pixel number of the map. 
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Fig. 2. — The upper limit of R, which corresponds to the instrument noise level, for the XPS 
to be able to successfully retrieve the CMB at the corresponding multipole number. Red, 
orange, green, light blue, blue and black color denote FWHM= 33, 24, 14, 10, 7, 5 arcmin, 
respectively, which then correspond to ESA Planck Surveyor LFI 30, 44, 70 GHz channel, 
HFI 100, 143 GHz, and the rest channels with higher frequency. The solid and dash lines 
are for LN = 1 and 100, respectively. 
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mon window function. The final power spectrum, after deconvolution, can be written as 
5T 2 = (2n)~ 1 £(£ + l)C e sm exp[b 2 £(£ + 1)], where C/ m is the (retrieved) CMB power spectrum 
and b = (v8 ln2) _1 FWHM. An error in deconvolution scale Ab produces an error in the 
power A(5T 2 ) = 2bAb£(£ + 1)5T 2 , so, for example, deconvolution with 1° FWHM on a 
smoothed map by 59' (i.e. overestimated by a mere 1 arcmin) shall overestimate the power 
at first Doppler peak by A(5T 2 2 20 )/5T 2 2 20 ~ 8.9% and the error goes with ~ £ 2 for higher 
£. Although in practice deconvolution involves estimating the inverse covariance matrix, it 
demonstrate that the window function is one of the most crucial issue in CMB power spec- 
trum estimation. The sensitivity of window function is also investigated in Sawangwit & 
Shanks (2010). 

Based on flat-sky approximation and Fourier transform, we can estimate the window 
function from bright point sources by taking two square patches of sky where there is a bright 
point source in the center: T ± (k) = ab(k) + c(k) + n(k), T 2 (k) = /3b(k) + c'(k) + n'(k), 
where k is the Fourier wavenumber, c, d and n, n' represent different CMB(+foreground) 
and noise, respectively. The CMB and noise are not correlated in different parts of the sky 
so they are now the "junk" in XPS whereas the bright point sources each in the center of 
the patches is manifestation of the inflight beam profile and thus are correlated signal. The 
window function is then = ^6 2 (fc) with a rescaling relation £ = 2nk/L, where L is the 
size of the patch. 

We simulate full-sky CMB signal with WMAP best-fit ACDM model on Planck 30 GHz 
channel and add 2 bright point sources on different parts of the sky, which is then convolved 
with beam 33 arcmin FWHM before adding noise with Ne = 10~ 2 /iK 2 . We extract 10° x 10° 
with the point source at the center. The point sources have amplitudes 40 and 46a, where a 2 
is the variance of the CMB. In Fig. 4 we present the simulation and the retrieving capability 
of the window function. One can see that XPS extended the estimate more than 10 dB. 



4. Window Functions of the WMAP Difference Assemly and Frequency Band 

Maps 

Having demonstrated that the window function can be retrieved with XPS from simu- 
lated bright point sources, we now employ the same method on WMAP DA maps: Kl, Kal, 
Ql, Q2, VI, V2, W1-W4. For each DA map, we extract the square patches where the top 
2 brightest point sources are in the center of the patche. 

Care has to be taken when one chooses the patch size. Fourier transform spreads the 
variance of the map into the power of different wavelengths. For stationary noise, choosing 
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Fig. 3. — Demonstration of the XPS retriving power and its limit. Top panel shows XPS 
(dots) of two maps with input spectrum Cg = 10 4 £ _2 (/zK 2 ) convolved with beam FWHM 
30 arcmin (dash-dot), but each added with different noise level (dotted): Ni = 1CT 1 and 
10~ 3 (/iK 2 ). Dash lines are the expected residual from XPS according to Eq.(6). Bottom 
panel shows binned (A£ = 10) power (plus sign) from average of 6 XPS of input WMAP 
best-fit ACDM convolved with W channel FWHM 12 arcmin (dash-dot), and added WMAP 
simulated noise (dotted). 
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Fig. 4. — Retrieving window function from bright point sources via XPS. Top panel shows 
the power spectra (solid red and blue curves) of 2 10° x 10° patches taken from a simulated 
full-sky map of Planck 30 GHz channel. In both patches CMB and a bright point source 
located at the center are convolved with 33 arcmin FWHM. Dash-dot black curve denotes 
the theoretical power spectrum of beam-convolved CMB from best-fit ACDM model whereas 
dash-dot red and blue the point sources with amplitude 47a and 40a, respectively, where 
a 2 is the variance of the CMB of the patch. Dotted line is the noise level 10 _2 /zK 2 and 
dash green line is the expected residual of the XPS from noise. In bottom panel we show 
the input window fucntion (solid curve) in dB and the retrieved window function (diamond 
sign), which is better than — 35dB as indicated in green line on top panel. 
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different patch area : changing from L\ to L 2 (with the same pixel size) shall shift the level 
of white spectrum by (Lq/L) 2 as the variance is fixed. For a point source, on the other 
hand, changing the patch size has an extra effect on the variance by (L /L) 2 , resulting in the 
large-scale power level of a point source being shifted by (L /L) 4 . Hence, choosing a smaller 
patch (L < Lq) enhances the relative power level of the point source to noise by (Lq/L) 2 , 
rendering a better profile of the window function. 

One cannot, however, choose as small sizes as possible. The first multipole number 
to display the profile of the window function is decided by the size of the patch L via 
£fi rs t = 27r/L, so one loses the large-scale profile if the patch size is too small. Furthermore, 
a fixed pixel size sets the maximum multipole number regardless of the patch size. So a 
smaller patch with a fixed pixel size has fewer Fourier modes with spacing A£ = 2n/L. 

We take the patch size L = 8, 5, 4, 3 and 3 degrees for K, Ka, Q, V and W band 
respectively and show the XPS results of all 10 DA maps in Fig 6: solid curves are the 
window functions from WMAP DA maps and diamond sign denotes the retrieved window 
functions from XPS. One can see in Fig. 5 that the XPS retrieves the window functions nicely, 
except for the W band 4 DA maps. The reason is that W band not only has the smallest 
beam size, but also higher noise level than the other bands. 

According to the previous XPS DAs results, our method for retrieving window function 
from bright point sources is demonstrated useful. We can now retrieve the window functions 
of the WMAP frequency band maps. The WMAP frequency band maps, together with the 
corresponding window functions, can be as useful as the DA maps in extracting the CMB 
power spectrum. The WMAP K, Ka, Q, V and W band have 1, 1, 2, 2 and 4 DA maps 
respectively, and the frequency band maps are produced from combining the DA maps, thus 
there is no direct measurement or estimation of the corresponding window functions for each 
map. One should note that the window functions of the DA maps at the same frequency 
band do not necessarily have the same profile, which is particularly true for W band (e.g. 
see Wl and W2 DA window functions in Fig. 5). 

In order to retrieve the window functions of the frequency band maps, we perform the 
same procedure as that on the DA maps in the previous section by choosing the same patches 
of the sky and patch size for the corresponding band. The retrieved window functions are 
showed in Fig. 6, and for comparison, we put the window function profiles of the DA maps 
for the corresponding band. Though there are 2 DA maps for Q and V band, the profiles 
are close to each other, and therefore one can see our retrieved window function for the 
frequency band maps are close to that of DA maps'. For W band, however, the difference 
in DA maps is pronounced in that the retrieved W band window function is compared with 
that of Wl DA map (solid curve) and W2 DA map (long dash). The profile of W3 and W4 
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are not plotted as they are close to that of W2 and Wl respectively. One can see that the 
window function profile of the W frequency band map is close to that of Wl (and W4), not 
W2 (and W3). 

5. Conclusion and Discussion 

In this paper we discuss the limit of the widely-used XPS and we demonstrate that it 
has the ability to decrease uncorrelated "junk" signal by a factor inversely proportional to 
the square root of the number of modes involved in the summation. It indeed can decrease 
the uncorrelated signal on small scales in the power spectrum, but has limited ability in 
eliminiating that of large scales. We then employ XPS to retrieve window functions via 
extragalactic point sources. One of the crucial steps for estimating CMB power spectrum 
is to know about the window function that convolves the observed signal. In this case the 
"junk" signal is composed of CMB and pixel noise. We take WMAP DAs as test samples 
first and employ XPS on extragalactic bright point sources to extract window function. We 
demonstrate it a useful method to extract window function without measuring planets for 
DAs. For frequency band maps combined from DA maps, we apply the same method to 
extract the window functions, which are unknown prior to our investigation. Our method 
will be useful for the upcoming Planck data. 

We acknowledge the use of Healpix 2 package (Gorski et al. 2005) and the use of 
Glesp 3 package. The author would like to thank Peter Coles and Dipak Munshi for useful 
discussions. 

REFERENCES 

Bennett, C.L. et al, 2003, ApJS, 148, 97 

Chiang, L.-Y., Naselsky, P.D., Coles, P., 2009, ApJ, 694, 339 

Doroshkevich, A.G. et al, 2005, IJMPD, 14, 275 

Gorski, K. M., Hivon, E., Banday, A. J.,Wandelt, B. D., Hansen, F. K., Reinecke, M., & 
Bartelmann, M., 2005, ApJ, 622, 759 



2 http : //healpix . jpl . nasa . gov/ 
3 http : //www . glesp . nbi . dk/ 



- 12 - 




Fig. 5. — The retrieved window functions of the 10 WMAP DA maps. Solid curves are the 
WMAP window functions for the DA maps and the diamond sign denotes the retrieved ones 
after crossing two patches of brightest point sources. One can see that the method is doing 
well on K, Ka, Q and V band DA maps, except for the W band DA maps, which is due to 
the small beam size and high level noise. 



- 13 - 





Fig. 6. — The retrieved window functions of the WMAP frequency band map. We employ 
the same procedure as that in retrieving the DA maps in Fig.5. The diamond sign denotes 
the retrieved window function and, for comparison, we put the window function profiles of 
the DA for comparison. W band is compared with Wl (solid curve) and W2 (long dashe) 
window function. 
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